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We analyze the behavior of supercooled fluids under shear both theoretically and numerically. 
Theoretically, we generalize the mode-coupling theory of supercooled fluids to systems under sta- 
tionary shear flow. Our starting point is the set of generalized fluctuating hydrodynamic equations 
with a convection term. A nonlinear integro-differential equation for the intermediate scattering 
function is constructed. This theory is applied to a two-dimensional colloidal suspension. The shear 
rate dependence of the intermediate scattering function and the shear viscosity is analyzed. We 
have also performed extensive numerical simulations of a two-dimensional binary liquid with soft- 
core interactions near, but above, the glass transition temperature. Both theoretical and numerical 
results show: (i) A drastic reduction of the structural relaxation time and the shear viscosity due 
to shear. Both the structural relaxation time and the viscosity decrease as A f~ v with an exponent 
v < 1, where j is the shear rate, (ii) Almost isotropic dynamics regardless of the strength of the 
anisotropic shear flow. 
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I. INTRODUCTION 

Complex fluids such as colloidal suspensions, polymer 
solutions and granular fluids exhibit very diverse rheo- 
logical behavior^], Shear thinning is among the most 
well-known phenomena. Such behavior is predicted for 
simple liquids as well0,0|, but the effect is too small to 
observe at temperatures well above the glass transition 
temperature. For supercooled liquids, however, the situ- 
ation is different. Recently, strong shear thinning behav- 
ior was observed by experiments performed on soda-lime 
silica glasses above the glass transition temper ature^. 
Yamamoto et al. have done extensive computer simu- 
lations of a binary liquid with a soft-core interaction^ 
near their glass-transition temperature and found non- 
Newtonian behavior. The same behavior was found for 
other systems such as Lennard- Jones mixtures 0, an d 
polymer meltsQ, too. For all cases, the structural relax- 
ation time and the shear viscosity decrease as A i~ v , where 
7 is the shear rate and v is an exponent which is less than 
but close to 1. For such systems driven far from equi- 
librium, the parameter 7 is not a small perturbation but 
plays a role similar to an intensive parameter which char- 
acterizes the "thermodynamic state" of the system fio). 
Such rheological behavior is interesting in its own right, 
but understanding the dynamics of supercooled liquids in 
a nonequilibrium state is more important because it has 
possibilities to shed light on an another typical and per- 
haps more important nonequilibrium problem, namely 
that of non-stationary aging. Aging is characterized by 
slow relaxation after a sudden quench of temperature be- 
low the glass transition temperature. In this case, the 
waiting time plays a similar role to (the inverse of) the 
shear rate. Aging behavior has been extensively studied 



in spin glasses (see Ref.^l an d references therein). Ag- 
ing is also observed in structural glasses There have 
been recent attempts to study the aging of structural 
glasses theoretically but no analysis and comparison 
to experiments or simulation results 12] have been pre- 
sented due to the complicated nonstationary nature of 
the problem. 

The relationship between aging and a driven, steady- 
state system was considered using a schematic model 
based on the exactly solvable p-spin spin glass by 
Berthier, Barrat and Kurchan|14j. This theory natu- 
rally gives rise to effective temperatures. The valid- 
ity of their idea was tested numerically for supercooled 
liquids 7J. There have also been attempts to observe ag- 
ing by exerting shear on the system instead of quenching 
the temperature [l^. Recently, there have been attempts 
to dev elop the mode-coupling theory for the sheared 
glasses [lfj Il7| . Fuchs and Cates have developed the 
mode-coupling theory for the sheared colloidal suspen- 
sions using projection operator techniques [17|. They 
have analyzed a closed equation for the correlation func- 
tion for a schematic model where the shear is exerted on 
all direction equally ("The Isotropically Sheared Hard 
Sphere Model"). In their model, shear is turned on at 
the initial time and therefore the dynamics are genuinely 
nonstationary. 

In this paper, we investigate the dynamics of super- 
cooled fluids under shear both theoretically and nu- 
merically for a realistic system. We extend the stan- 
dard mode-coupling theory (MCT) for supercooled flu- 
ids and compare the solutions with the numerical simu- 
lation results. We mainly focus on the microscopic origin 
of the rheological behavior. The relationship with the 
more generic aging problem will not be discussed here. 
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Since our goal is the investigation of shear thinning be- 
havior, we have neglected violations of the fluctuation- 
dissipation theorem. We start with generalized fluctu- 
ating hydrodynamic equations with a convection term. 
Using several approximations, we obtain a closed non- 
linear equation for the intermediate scattering function 
for the sheared system. The theory is applicable to both 
normal liquids and colloidal suspensions in the absence 
of hydrodynamic interactions. Similar approach has been 
applied to the self-diffusion of the hard sphere colloidal 
suspension at relatively low densities by Indrani et al. |l8| . 
Extensive computer simulations are implemented for a 
two-dimensional binary liquid interacting with soft-core 
interactions. The effects of shear on microscopic struc- 
ture, the structural relaxation time, and rheological be- 
havior are discussed. Results both from the theory and 
the simulation show good qualitative agreement despite 
of the differences between the systems considered. Spe- 
cial attention is paid to the directional dependence of the 
structural relaxation. 

The paper is organized as follows: In the next section, 
we develop the MCT for sheared suspensions. Complex- 
ities which do not exist in mean-field spin glass models 
and how those complexities should be treated are elu- 
cidated here. A possible way to explore the situation 
without invoking the fluctuation-dissipation theorem is 
also discussed. In Section III, the model and our simula- 
tion method are explained. The results both from theory 
and simulation are discussion in Section IV. In Section V 
we conclude. 



interactions 22] but we shall neglect these effects. The 
friction term is specific for the colloidal case. In the case 
of liquids, it should be replaced by a stress term which 
is proportional to the gradient of the velocity field mul- 
tiplied by the shear viscosity. Both cases, however, lead 
to the same dynamical behavior on the long time scales 
which are of interest here. f#(r,i) is the random force 
which satisfies the fluctuation-dissipation theorem of the 
second kind (2nd FDT) 23]; 

(/fl,,(r,<)/fl J (r',i'))o = 2fc B Tp(r,i)CoM(r-r')<5(i-i') 

. (IL3) 

for t >t' , where (■ ■ • )o is an average over the conditional 
probability for a fixed value of p(r, t) at t = t'. Note that 
the random force depends on the density and thus the 
noise is multiplicative. This fact makes a mode-coupling 
analysis more involved as we discuss later in this section. 
We assume that the 2nd FDT holds even in nonequilib- 
rium state since the correlation of the random forces are 
short-ranged and short-lived, and thus the effect of the 
shear is expected to be negligible. The first term on the 
right hand side of the equation for the momentum is the 
osmotic pressure term. Here we assume that the free en- 
ergy T is well approximated by that of the equilibrium 
form and is given by the well-known expression 



(II.4) 



0F~ /dr p(r){lnp(r)/ Po -l} 

dr t /dr 2 6p(r 1 )c(r 12 )6p(r 2 ) 7 



II. MODE-COUPLING THEORY 

We shall consider a two-dimensional colloidal suspen- 
sion under a stationary shear flow given by 



v (r)=r.r= (jy,0), 



(11.1) 



where 7 is the shear rate and (T) a p = jS ax Si3 y is the ve- 
locity gradient matrix. Generalization to higher dimen- 
sions is trivial. We start with the generalized fluctuating 
hydrodynamic eauations |19t l20| . This is a natural gener- 
alization of the fluctuating hydrodynamics developed by 
Landau and Lifshitzj2lJ to short wavelengths where the 
intermolecular correlations becomes important. Fluctu- 
ating fields for the number density p(r, t) and the velocity 
v(r, t) of the colloidal suspension obey the following set 
of nonlinear Langevin equations 



dt 



TTlV 



8J~ 

H = -P^s Cop(v - Vq) + ffi, 

dp 



(II.2) 

where m is the mass of a single colloidal particle, T is 
the free energy of the system and C$ is a bare collective 
friction coefficient for colloidal particles. Co has a weak 
density and distance dependence due to hydrodynamic 



where — l/k B T, pg is the average density, and c(r) is 
the direct correlation function. We have neglected cor- 
relations of more than three points, such as the triplet 
correlation function 03(1*1, r 2 , ra), whose effect becomes 
important for the fluids with stronger directional inter- 
actions such as silica[24j. Under shear, it is expected 
that c(r) will be distorted and should be replaced by 
a nonequilibrium steady state form c NE (r), which is an 
anisotropic function of r. It is, however, natural to ex- 
pect that this distortion is small on the molecular length 
scales which play the most important role in the slowing 
down of structural relaxation near the glass transition. 
The distortion of the structure under shear is given up 
to linear order in the shear rate bv|25j. 



S NE (k) = S(k) { 1 



k • r ■ k dsw 1 



2fcD n dfc 



J 



(II.5) 



where S(k) and 5 NE (k) is the static structure factor in 
the absence and in the presence of the shear, respec- 
tively. Do = /sbIYCo is the diffusion coefficient in the 
dilute limit, k = |k|, and k = k/|k|. The direct cor- 
relation function in the Fourier representation of c(r) is 
related to the structure factor by nc(fc) = 1 — 1/S(k). 
From ea. HII.5| l. we find that the distortion due to the 
shear is characterized by the Peclet number defined by 
Pe= ja 2 /Do, where a is the diameter of the particle. 
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Hereafter we shall neglect the distortion and use the di- 
rect correlation function at equilibrium, assuming the 
Peclet number is very small. 

We linearize ea. l|II.2(l around the stationary state as 
P — Po + &P an d v = v + Sv, where po is the average 
density. Transforming to wave vector space, we obtain 
the following equations; 



d_ 

at 



k r 



d_ 
dk 



dt dk 



ik 



pS(k) 



Spk(t) 



Spk(t) = — J k (t), 
m 



k r 

1 

m(3 



kj Mt) 



--Mt) + fm(t), 
m 

(II.6) 

where Jk(i) = mpok ■ <5vk(i) is the longitudinal mo- 
mentum fluctuation, and J = J dq/(27r) 2 . We have 
neglected the quadratic terms proportional to J q Jk-q- 
Note that ea. Hll.6l) does not contain coupling to trans- 
verse momentum fluctuations even in the presence of 
shear. 

The direct numerical integration of eq. (jll.6(l is in prin- 
ciple possible but it is expensive and not theoretically 
enlightening 26] . Rather we shall construct the approx- 
imated closure for the correlation functions, a so-called 
mode-coupling approximation!^ |2||,|2il. There are sev- 
eral approaches to derive mode-coupling equations, in- 
cluding the use of the Mori-Zwanzig projection operator 
and a decoupling approximation 27], or implementation 
of a loop expansion developed in the context of the equi- 
librium critical phenomena and generalized to the dy- 
namic case[29l l30j. Both approaches lead to essentially 
the same equations if the system is at equilibrium. Un- 
der noncquilibrium conditions, however, the loop expan- 
sion approach is more straightforward and flexible. We 
shall adopt the loop expansion approach to the nonlinear 
Langevin equation with both multiplicative noise and the 
full convection term. Ea. l|II.6|l can be cast to a general 
form of the nonlinear Langevin equation written as 



~dt 



-zVijkXjXk + Jr.. 



(11.7) 



where x(t) = (5pk(t), Jk(t)) is a field variable, pij is the 
linear coefficient matrix, and the nonlinear coupling co- 
efficient Vijk is the vertex tensor which satisfies the sym- 
metric relation Vy& = Vikj- Finally, fii,i(t) is the random 
force field which satisfies the 2nd FDT, 



(/R,i(t)/«j(«')>o - k B L rj {*)5{t - t'), 



(II. 



where Ly(x) is the x-dependent Onsager coefficient 
which is to be expanded to the lowest order as 



L ij (x)=Lf+L ( V kXk 



(II.9) 



where L 



(i) 

ij.k 



dLij(x)/dxk\x=o- Comparing ea. HII.7|) 
with eq. 1)11.6(1 , the elements of the linear coefficient ma- 
trix are given by 



/VkPk' 
PJ k J^ 



d 

k r • — <s kik , 

ok 

ik 

— Ok.k' 



ik 



f3S{k) 



'k,k' 



k • r — - k • r • k) 5 k .k< - -<5 k ,k< 

ok I m 



(11.10) 

Nonzero elements of the vertex tensor are given by 



V 



-fkPk'Pk" 



■^{ik.k'c(fc')+*k-k"c(fc")} 



^k,k'+k", 
(11.11) 

(11.12) 



where V is the volume of the system. Finally, 
L JkJ k ' = CoPoTVSk-k' 

L JuJw,P k » = T Co<W,k". 
All other components are zero. In the above expressions, 
_ (2^) 2 



Jk,k' 



V 



-S(k-k') 



is the Dirac delta function. We construct the clo- 
sure equation for the correlation function Cij(t,t') = 
(xi(t)xj(t')) . Since we are treating the stationary 
state, the time translation invariance holds and, thus, 
dj(t,t') = Cij(t — t'). A general scheme for the loop 
expansion method for the Langevin equation with mul- 
tiplicative noise has been discussed by Phythian[3lJ. Up 
to one-loop, the equation for the correlation function is 
written as 32] 

= / dh E iQ (t-ti)C Q j(ti -t') 

J —oc 

+ f Ah D ia (t - t^Gifi! -t') 
dG tJ (t-t>) 



dt 

t 



p ia G aj (t - t') = 5 i3 5{t - t 1 ) 



dt\ SjQ,(i — ti)G a j(ti — t) 

with the memory kernels defined by 

E.y(i) =Vi a pG a \(t)Cp IJ/ (t)V\ l ij 

+ kBV ia0 G a x(t)G Pll (t)L l Vj 

Dij(t) =-^Vial3G a \(t)C/3 f _ l (t)Vj\ fl 

+ 2k B V iaf3 G aX {t)C 0ll {t)Lfl^ 
+ Zk B L%& aX {t)C p ,{t)V 3 ^. 



(11.13) 



(11.14) 
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In these expression, we have introduced the propagator 
defined by 



»(*" 



Gji(t'- 



Sxj (t) 



(11.15) 



-t) is the conjugate of the propagator. 
Ea. (|II.13|) together with ea. (jIL14|) are so-called mode- 
coupling equations. 

Without further assumption, one needs to solve the 
coupled equation for the correlation function and the 
propagator. The propagator can be eliminated, however, 
if the fluctuation-dissipation theorem of the first kind 
(1st FDT) is valid. The 1st FDT relates the correlation 
function to the response function, \ij (t) via the equation 



Xij (t) 



0{t) dC l3 (t) 
knT At 



(11.16) 



where 6{t) is the heaviside function. Note that the re- 
sponse function Xij(t) is generally neither the same as, 
nor proportional to, the propagator Gij(t). Xij(t) repre- 
sents the response of the system to the perturbation via 
external fields or through the boundary, whereas Gij(t) 
represents the response to a random force field. If the 
Langevin equation is linear, both are the same, but this 
is not true in general for the nonlinear Langevin equation. 
To one-loop, to be consistent with the mode-coupling ap- 
proximation, the response function is written as|32l l33j 
as 



Xi j (t-t') = ^G ik (t-t'){M^ + L^} 



+ 7^ J&tl Gik(t - ti)Vka/3G a \(ti — t')C 0fJ ,(tl - t') 

x{M«+Z«} A:WI , 

(11.17) 

where the tensors Mjp and are defined by 

Mij(x) ~ M[f + M^\x k . Mij(x) is the matrix which is 
defined by the reversible part of the nonlinear Langevin 
equation. In general, we can express the Langevin equa- 
tion as 



~df 



ss ss 

r,(xS +M ij {x)— +L ij {x)— (11.18) 



where 5* is the entropy of the whole system and Vi (x) is a 
term which originates from the nonequilibrium constraint 
such as the convection for the sheared system. From 
this definition, the reversible term does not contribute 
to the entropy production and thus the matrix M i3 (x) 
is antisymmetric. For the system considered here, the 
non-zero elements are given by 



M, 



PkJ k 



= -Mj • , Pk = -ikT (p Q V5 k - k > +<5 / o k +k') 



(11.19) 

If the system is in equilibrium, one may eliminate 
the propagator in favour of the correlation function in 



ca. llTLh!) by using the 1st FDT, ea. fjTT6l) |3jJ |3J, |35j . 

If the system is out of equilibrium, one has to solve the 
coupled equations 1)11. 13|1 . For the case of the p-spin- 
glass model with an external drive, Berthier et al. have 

They have 



analyzed the equation similar to l|II.13l) 14]. 
found that there is a systematic deviation from the 1st 
FDT which reminiscent of violations of FDT that occur 
during aging. Extensive computer simulation supports 
these results 7J. A similar analysis for real fluids is nec- 
essary but it is inevitably more involved due to the com- 
plicated tensorial nature of the non-linear coupling. In 
this paper, we shall not focus on the fundamental prob- 
lem of the validity of the fluctuation-dissipation theorem 
and assume that the 1st FDT is valid. We shall show 
that even with this simplification, the theory can explain 
qualitative several aspects of the dynamical behavior very 
well. Using the 1st FDT, as in equilibrium case, one can 
eliminate Gg( t) from ea. (|II.13|) and the final expression 
is given by 32] 



dh M ik {t-h)C kj {ti-t') 



with the memory kernel given by 



(11.20) 



Mijit) = -\v ia0 C aX {t)C 0ll {i) {V feV - 2£ fcV } £7,7/(0). 

(11.21) 

These equations are almost identical to the conventional 
mode-coupling equations |34l 1351 ] except for (—2Ckx^) 
term which appears in the vertex term of the memory 
kernel. The third term on the left hand side of ea. l|IL20|) 
enters naturally to guarantee the time-reversal symmetry 
of the correlation function, djk is defined by 



r _r(0)o(3) , r(l) o(2) , r (1) o( 2 ) CTT 99^1 



where we have defined S^' and by 



_ 

13 dxidxj 



and S 



x=0 



(3) 



d 3 S 



dxidxjdxk 

x=0 

(11.23) 

The term (— 2£/ c- \ Ai ) is a new term which originates from 
the multiplicative noise. Since the density- and velocity- 
dependent part of the entropy is given by 



— H T 

2Tp{v) T 



(11.24) 
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and S^l are given by 



? (2) 



0(2) 
w PkP k / 



C(3) 
PkPk'Pk" 



o(3) 

D J k J k ,p k „ 



1 



mpoTV 
k B 



Sk,-k> 



PoS(k)V 

~ 2 T/9 Ok+k'+k",0 



(11.25) 



1 



?k+k'+k",0- 



Combining this with ea. (|II.12fl . we find that Cijk = for 
all elements. Therefore, the final expression for the mode- 
coupling equation in the present case becomes equivalent 
with the conventional onej^- Note that the situation 
will change if we consider different physical situation: For 
example, if we start with the diffusion equation, which is 
obtained in the overdamped limit of ea. 1)11.6(1 . Cijk plays 
an essential role|3^. 

Now let us substitute all matrix elements given by 
cas. (|11.10[) and (|ll.ll|l to ca. (|11.20() . we have the set of 
equations given by 



at dki cpkp "' 

d_ 

at 



(t) = -ikC JkPk ,(t) 



k r 

ik 

'mpS(k) 
ft 



d_ 
dk 



k-r-k) c JkPk ,(t) 
Co, 



CpkP k ' (t) - —Cj kPk , (t) 



(11.26) 



- / dh V<5C(k,k",t-t 1 )C/ k „p k ,(ii) 
71 J o t „ 



where 



E 



V 



(27T) 2 



dk. 



(11.27) 



The memory kernel <5£(k, k', t) is given by 
S((k,k',t) 

EEEEE 

^kPqiPqa Cp qi P P1 

qi q2 pi P2 k" 
x^ k „ Ppi p P2 ^, Jk ,(0) 



m 
~2 



2 PP2 



(*) 



2nV 



EE^ PqPk-q^PqPp WCpk-qP-k-pW 



q p 

x V. 



f 



^-k'P P P-k'- P 

^ k (q,k-q)C PqP _ p (<)C Pk _ qP _ k , +p (i) 



2n*Nl3 7 q 7 P 

xVk'{p,k' -p), 

(11.28) 

where N is the total number of the particles in the system 
and we have defined the vertex function 

V k (q,k-q)=k-qnc(g)+k-(k-q)nc(|k-q|). (11.29) 





FIG. 1: (a) Geometry of shear flow, (b) Shear advection in 
real space, (c) Shear advection in Fourier space. 



In the derivation of eq. 2Rf> , we have used the momen- 
tum equal-time correlation function given by 



mnV 



(11.30) 



If the system is in equilibrium, time and space transla- 
tional invariance is satisfied and the correlation function 
becomes C a Mp( T ') (t) — Ca(r-r')/3(o)(*) or ; m k space, 
Ca k /3 k , (t) = Ca k /3_ k (t) x (5k, -k'- In the presence of shear, 
however, this should be modified as[37j. 

C Q (r)/3(r')(*) = Ca(r-r'(t)),3(0)( i ) I (11.31) 

where we defined the time-dependent position vector by 

r{t)=exp[Tt]-r = r + jtye x , (11.32) 

where e x is an unit vector oriented along the x-axis. In 
wave vector space, this is expressed as 



(11.33) 



C a k f3 k ,(t) - C ak(t) fj_ k {t) X <5 k(t ) _ k / 
= CW'M/VM X ^k-k'(-t) 

with the time-dependent wave vector defined by 

k(t) = exp['lY] • k = k + jtk x e y , (11.34) 

where 'T denotes the transpose of I\ Figure 1 (b) 
shows how the shear flow advects a positional vector 
r by ea. l|II.32l) in a time interval of duration t. The 
corresponding time-dependent wave vector, ea. (|II.34I . is 
shown in Fig.l (c). Eas. i|II.31() and I|II.33J| state that the 
fluctuations satisfy translational invariance in a reference 
frame flowing with the shear contours. Using this prop- 
erty, ea. ljll.28(l becomes 

<5C(k,k',i) 
1 



x Vt( t )(p(*),k(i) - q(t))$ p ,q(t)<Jk',k(t) 

2plNV(3 j k ^ q ) C 'p< 5 P-q( t ) ( i ) C 'pk-qP-k( t ) + q( t ) M 

x Vk(t)(p(t),k(t) - q(t))5 k /,k(i) 
5((k,i) x 5 WMt) . 

(11.35) 
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Introducing the intermediate scattering function by 
F(k,t) ee ± Cpk( _ t)P _ k (t) = hs Pk ^ t) (t)5p. k (0)), 



(11.36) 



5((k, t) can be rewritten as 

5C(k,t) = -L [ V k (q,k-q)F(q(t),t)F{k(t)-q(t),t) 

x v m { P {t)Mt)-<i(t)). 

(11.37) 

This is the memory kernel in the presence of the shear. 
This has exactly the same structure as the equilibrium 
one|27j except for the time dependence appearing on the 
wave vectors. The physical interpretation of ea. (|II.37|) 
is simple: The memory kernel has the general structure 
of two correlation functions sandwiched by two vertex 
functions. This means that the interactions of two fluc- 
tuations with modes q and k q scatter at a certain time 
and then propagate freely in a "mean field" and after a 
time t, they recollide and interact. Under shear, how- 
ever, by the time the second interactions take place, the 
fluctuations are streamed away by the flow. 
We also define the cross correlation function 



H(k,t) ee ± Cjk( _ t)P _ k (t) = l<j k( _ t) (t)apj(o)>. 



(11.38) 



Then, ea. (|II.26|> can be rewritten as 



dF(k, t) 
dt 

dl 



= -ik(-t)H(k,t). 

= -k(-t) •r-k(-t)ff(k,t) 



(11.39) 



mPS{k(-t)) v ' ' 7n 

1 r* 

— / dt' 6((k(-t),t-t')H(k,t'). 
m Jo 



Note that in the above equation, the differential operator 
k r • d/dk disappears because 



dF(k,<) dF(k,t) 



dt 



dt 



■k(-t)-r. 
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dk(-t) 



F(k,t). (11.40) 



The memory kernel eq. (|II.37|) appearing in eq. (|II.39|) con- 
tains the nonlinear coupling with the correlation func- 
tion itself and therefore the equation should be solved 
self-consistently. Indrani and Ramaswamy have derived 
an equation similar to ea. ljll.39j) for the velocity corre- 
lations of a single tagged particle in a three-dimensional 
hard sphere colloidal suspension |Tsj. They were inter- 
ested in the relatively low density regime, and did not 
solve the resulting equation self-consistently. 

For colloidal suspensions the relaxation time of the mo- 
mentum fluctuations is of the order of r m = to/Co an d is 
much shorter than the relaxation time for density fluctua- 
tions which is of the order of or longer than Td = a 2 /Do. 
In other words, we may invoke the overdamped limit. 



Thus, we neglect the inertial term dH(k, t)/dt. Likewise 
the first term of the right hand side in the second equa- 
ls estimated to be of order of jT„ 



tion of ea. l|II.39|l is estimated to be of order of 7T m and 
thus should be very small as long as the Peclet number 
is small. Thus, the equation for the momentum fluctua- 
tions may be written as 



= - SfTT ^ HK t) - ^H(k, t) 
m(3S{k(-t)) y ' ' m K ' ; 

1 f* 

/ dt' 5C(k(-t),t-t')H(k,i/). 

m Jo 



(11.41) 



Substituting this back into the first equation of (|II.39|) . 
we arrive at the equation for the intermediate scattering 
function: 



dF(k,i) _ Dok(-tf 



dt 



F(k,t) 



S(k(-t)) 
df'M(k(^),t-t')^M : 



(11.42) 



where the memory kernel is given by 



M(k, t) 



Do k 



2 P o k(t) Jq 
xF(k(t)-q(t),i)F(q(t),i) 



Vt(q,k-q)V k(t) (q(t),k(t)-q(t)) 



(11.43) 

are the final mode-coupling ex- 



Eas. pL42fl and pL43|l 

pressions. We have started from the equation for both 
the density and the momentum fields and proceeded via 
the standard mode-coupling approach, taking the over- 
damped limit at the end. It should be possible to de- 
rive the same result starting from the diffusion equation 
with an interaction term which is obtained by taking the 
overdamped limit in ea. l|II.6|) . In order to arrive at the 
same result, however, one needs to use a different resum- 
mation scheme that involves the irreducible projection 
operator^. 

Eq. l|II.42|) is a non-linear integro-differential equation 
that can be solved numerically. An efficient numerical 
routine to solve the mode-coupling equation is elucidated 
in Ref . |39l| . Since our equation is not isotropic in k due 
to the presence of shear, the numerics are more involved 
than in the equilibrium case. We shall solve the equation 
by dividing the lower half plane into an Nk x (Nk —l)/2 
grid, where Nk is the grid number which we have cho- 
sen to be an odd number. We do not need to consider 
the upper half plane because it is a mirror image of the 
lower one. k x , for example, is discretized as k x .o — 
-k c ,k x ,i = -k c + 5k,- ■ ■ ,k x ,N k = -k c + N k 5 k = k c , 
where 8k = 2k c /Nk is the grid size and k c is a cut-off wave 
vector. Any value outside the boundary, \k x \, \k y \ > k c , 
is replaced by the value at the boundary. 



III. SIMULATION METHOD 

To prevent crystallization and obtain stable amor- 
phous states via molecular dynamics (MD) simulations, 
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we choose a model two-dimensional system composed of 
two different particle species 1 and 2, which interact via 
the soft-core potential 

v ab (r) = eK fc /r) 12 (III.l) 

with cr a b = (a a + f&)/2, where r is the distance between 
two particles, and a, b denote particle species (s 1, 2). We 
take the mass ratio to be m 2 jm,\ = 2, the size ratio to be 
02/ci = 1-4, and the number of particles TV = Ni + N 2 , 
Ni = N2 = 5000. Simulations are performed in the 
presence and absence of shear flow keeping the particle 
density and the temperature fixed at n = ri\ + n 2 — 
0.8/o-f (ni = N x /V, n 2 = N 2 /V) and k B T = 0.526e, 
respectively. Space and time are measured in units of G\ 
and to = (miCTj/e) 1 / 2 . The size of the unit cell is L = 
118eri. In the absence of shear, we impose microcanonical 
conditions and integrate Newton's equations of motion 

at m a at 

after very long equilibration periods so that no apprecia- 
ble aging (slow equilibration) effect is detected in various 
quantities such as the pressure or in various time corre- 
lation functions. Here, r? = [r a xl ,r a yl ) and p? = (p xi ,Py t ) 
denote the position and the momentum of the «-th parti- 
cle of the species a, and ff is the force acting on the «-th 
particle of species a. In the presence of shear, by defin- 
ing the momentum p'° = p° — m a ^T ui : e. x (the momentum 
deviations relative to mean Couette flow), and using the 
Lee-Edwards boundary condition, we integrate the so- 
called SLLOD equations of motion so that the tempera- 
ture ksT (= W^^aEitP'n 2 /^) is kept at a desired 
value using a Gaussian constraint thermostat to elimi- 
nate viscous heating effects flo). The system remains at 
rest for t < for a long equilibration time and is then 
sheared for t > 0. Data for analysis has been taken and 
accumulated in steady states which can be realized after 
transient waiting periods. 

We shall calculate the incoherent and the coherent 
parts of the scattering function for the binary mixture 
by using the definitions 

F,(k, t ) = _L/ V el-W-^W-^^A (III.3) 

A a \ , / 
' i— 1 ' 

and 

. , N a N h . 

F ah {k,t) = V^ e H k H)^W]^ e N(»)]\ 
\ i=i j=i I 

(III.4) 

respectively with a,b € 1,2. The a relaxation time r Q of 
the present mixture, which is defined by 

F n (k ,T Q ) ~ F s (k 0) r a ) = e" 1 , (III.5) 

is equal to T a ~ 1800 time units in the quiescent state 
for |ko| = 2n/ai. 



IV. RESULTS 
A. Microscopic Structure 

The partial static structure factors S , a h(k) are defined 

as 

S ab (k) - J dr e 4k - r (n a (r)n b (0)), (IV.l) 

where 

N a 

fi„(r) = ^(r-rj) (a = 1,2) (IV.2) 
i 

is the local number density of the species a. Note the 
dimensionless wave vector k is measured in units of erf . 
For a binary mixture there are three combinations of par- 
tial structure factors, 5n(k), ^22(k), and S^k). They 
are plotted in Fig. 2 (a) in the quiescent state after taking 
an angular average over k. A density variable represent- 
ing the degree of particle packing, corresponding to the 
density of an effective one component system, can be de- 
fined for the present binary system by 

p cS (r) = ^nx(r) + o\h 2 {v). (IV.3) 

The corresponding dimensionless structure factor is given 
by 

S pp (k) = af 4 J dr e 4k - r (,5p cff (r)^ cff (0)) 

= n 2 S n (k) + n 2 2 (a 2 /a 1 ) i S 22 (k) + 2n 1 n 2 (a 2 /<r 1 ) 2 S 12 (k), 

(IV.4) 

where 5p c g — p c s — (p c s). One can see from Fig. 2 (b) 
that S pp (k) has a pronounced peak at k ~ 5.8 and be- 
comes very small (~ 0.01) at smaller k demonstrating 
that our system is highly incompressible at long wave- 
lengths. Because S pp (k) behaves quite similarly to S(k) 
of one component systems, we examine space-time corre- 
lations in /5 c ff(r) rather than those in the partial number 
density n a (r) for the present binary system. The usage 
of Poff(r) makes comparisons of our simulation data with 
the mode-coupling theory developed for a one component 
system more meaningful. 

We next examine the anisotropy in the static structure 
factor S pp (k) in the presence of shear flow. Figures 3-6 
show 5*p P (k) plotted on a two-dimensional k x — k y plane 
(upper part) and the angular averaged curves (lower 
part) within the regions (a)-(d) obtained at 7 = 10~ 4 , 
10~ 3 , 10 -2 , and 10 _1 , respectively (in units of t^ 1 ). One 
sees that, at the lowest shear rate (7 = 10~ 4 ), the shear 
distortion is negligible but at higher shear rate the dis- 
tortion becomes prominent. At 7 = 10 -3 , at all regions 
except for the region (d), the peak heights of S pp (k) start 
decreasing. This is contrary to the estimate from the lin- 
ear response theory given by ea. ljll.5(l which predicts no 
distortion in the region (a) and (c) but distortion to the 
higher, (b), and lower, (d), peaks. This seems to indicate 








(b) 






k=10 






1 


k=2.9 













10 



FIG. 2: Partial structure factors S a b(k) in (a) and S pp (k) in 
(b) denned by ea. ljIV.4(l for the present binary mixture. 



that the nonlinear effects due to the shear become impor- 
tant. Ronis has explored the higher shear regime for the 
structure factor of hard sphere colloidal suspensions and 
concluded that at higher shear, the peak should be al- 
ways lower than the equilibrium value together with a 
shift that depends on the direction plf. Figure 6 shows 
that peaks in all directions have been lowered and the 
peak with maximal distortion (region (b)) is shifted to 
the lower wave vectors while the opposite is true for the 
shift of the region with minimal distortion (region (d)). 
The qualitative agreement with Ronis' theory is good 
but it is not clear that our results can be explained by a 
simple two-body theory such as that of Ronis. Recently 
Szamel has analyzed S(k) for hard sphere colloidal sus- 
pensions up to the linear order in 7^2- He took the 
three-body correlations into account and found quantita- 
tive agreement with the shear viscosity evaluated using 
S(k). 

It should be noted that because our system is a liquid, 
we cannot directly refer to the Peclet number since the 
bare diffusion coefficient Do does not exist. Therefore, 
a direct and quantitative comparison with the theories 
discussed above is not possible. However, estimates from 
the relaxation times self-diffusion coefficient evaluated in 
Ref. @ allow us to estimate a Peclet number in the range 
between 10 _1 and 10 2 , which corresponds to the highest 
shear rates explored in the theoretical analysis of this 
paper. 



B. Intermediate Scattering Function: Numerical 
Results 



Here, we examine the dynamics of the local density 
variable /5 e g(r, t). To this end, we defined the intermedi- 
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FIG. 3: S pp (k) for 7 = 10~ 4 . The solid line is S pp (k) at 
equilibrium. Dots represents those observed in the region in- 
dicated in the (k x , k y ) plane above. 



ate scattering function 
F pp (k,t) =n 2 Fn(k,t) 



n 2 (a 2 /CTi) 4 F 22 (k,i) 



(IV.5) 



+ 2nm 2 {(T2/ai) 2 F 1 2(k,t) 

by taking a linear combination of the partial scattering 
functions defined in ea. (|III.4| . Note that F pp (k, 0) = 
5 pp (k) by definition. To investigate anisotropy in the 
scattering function F pp (k, t) , the wave vector k is taken 
in four different directions kio, kn, koi, and k_n, where 



(iv.6) 



and /x, v G 0,1 as shown in Fig. 7. The wave vector k 
(in reduced units) is taken to be 2.9, 5.8, and 10 (see 
also Fig. 2 (b)) Because we use the Lee-Edwards periodic 
boundary condition, the available wave vectors in our 
simulations should be given by 



2tt 



(ne x , (m - nD x )e y ), 



(IV.7) 



where n and m are integers, D x = L-yt is the differ- 
ence in x-coordinate between the top and bottom cells 
as depicted in Fig. 6. 5 of Ref. [40t| . To suppress statistical 
errors, we sample about 80 available wave vectors around 
k, i) using eas. ljIV.5|) and l|III.4|l . 



and calculate F pp 
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FIG. 4: S'pp(fc) for 7 = 1CT 3 . All symbols are as in FigOU 



FIG. 6: S pp (k) for 7 = 1CT 1 . All symbols are as in Figg] 



and then we average F pp (k,t) over sampled wave vec- 
tors. The sampled wave vectors are shown as the spots 
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FIG. 5: S pp (k) for 7 = 1CT 2 . All symbols are as in FigEU 



FIG. 7: Sampled wave vectors. 



in Fig. 7. 

Figures 8, 9, and 10 show F pp (k, t)/S pp (k) for k = 2.9, 
5.8, and 10, respectively. We have observed that all of 
the partial scattering functions Fu(k, t), Fi 2 (k, i), and 
F22(k, t) behave in a similar manner as F pp {k, t) 1 which 
demonstrates that the effective single component scatter- 
ing function typifies the dynamics of the whole system. 
Several features are noticeable. First, the quantitative 
trends as a function of k are similar for different values 
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of 7. Secondly, shear drastically accelerates microscopic 
structural relaxation in the supercooled state. The struc- 
tural relaxation time r Q decreases strongly with increas- 
ing shear rate as r a ~ A f^ v with v ~ 1. Lastly, the 
acceleration in the dynamics due to shear occurs almost 
isotropically. We observed surprisingly small anisotropy 
in the scattering functions even under extremely strong 
shear, 7r Q ~ 10 3 . A similar isotropy in the tagged par- 
ticle motions has already been reported in Ref. 6]. The 
observed isotropy is more surprising than that observed 
in single particle quantities. In particular, the fact that 
different particles labels are correlated in the collective 
quantity defined in ea. (|III.4fl means that a simple trans- 
formation to a frame moving with the shear flow can- 
not completely removed the directional character of the 
shear. Our results provide post facto justification for the 
isotropic approximation of Ref. 17]]. This simplicity in 
the dynamics is quite different from behavior of other 
complex fluids such as critical fluids or polymers, where 
the dynamics become noticeably anisotropic in the pres- 
ence of shear flow. 




FIG. 8: F (k, t)/S(k) at ko\ = 2.8 for various shear rates and 
at the different observing points (a) kio, (b) kn, (c) koi, and 
(d) k_n as explained in Fig0 



C. Intermediate Scattering Function: MCT Results 

We evaluate F(k, t) for the two-dimensional colloidal 
suspension theoretically using ea. l|II.42|l with ea. l|II.43|l . 
Following the procedure explained in Section [HJ we have 
solved the equations ea. (|II.42|) self-consistently. For the 
static correlation function c{k) and S(k), the analytic ex- 
pressions derived by Baus et al. were used (see Appendix 
A) |43| . The number of grid points was chosen to be 
Nk = 55. We have also calculated solution of eo. l|II.42(l 
for the larger grid sizes up to N k — 101 but qualita- 
tive differences between results obtained with different 




FIG. 9: F(k,t)/S(k) at kai = 5.8. All symbols are as in 
FigEl 




t 



FIG. 10: F(k,t)/5(k) at ken = 10. All symbols are as in 
FigEJ 



grid numbers were not noticeable. We chose k c o~ = 10n 
as the cut-off wave vector. The best estimate for the 
transition density is tfi c — 0.72574 for Nk — 670 at equi- 
librium, where <fi = 7r<7 2 po/4 is the volume fraction. For 
Nk = 55, a higher value, </> c = 0.76645 is obtained. Fig- 
ures 11-13 show the behavior of -F(k, i) renormalized by 
its initial value S(k) for <f> — <fi c x (1 — 10~ 4 ) for var- 
ious shear rates from Pe=10~ 10 to 10 _1 , where Pe is 
the Peclet number defined by Pe= 7cr 2 /-Do- The wave 
vectors were chosen to be ka = 3.0 (Fig. 11), ka = 6.7 
(Fig. 12), and ka = 12.0 (Fig. 13). ka = 6.7 is close to 
the position of the first peak of S(k). At each wave vec- 
tor, we have observed F(k,t)/ S(k) for four directions 
denoted by (a)-(d) in Fig. 7 and defined by ea. l|IV.6|) . For 
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FIG. 11: F(k,t)/S(k) for ka = 3 for different observing 
points and for various shear rates. From the right to the left, 
Pe = 10" 10 , 10" 7 , 10" 5 , 10" 3 , and 10" 1 . The thick dotted 
lines are at (a) in Fig|7] The thin dotted lines are at (b). The 
thick solid lines are at (c). The thin solid lines are at (d). 
The density is (j> = (f> c x (1 — 1CF 4 ). The time t is scaled by 
a 2 /D . 




10" 4 10" 2 10° 10 2 10 4 10 6 10 8 10 10 10 12 
t 



FIG. 12: F(k,t)/S(k) for ka = 6.7. All symbols are as in 
Figim 



shear rates smaller than Pe = 10~ 10 no effect of shear is 
observed. For Pe > 10~ 10 , we observe a large reduction 
of relaxation times due to shear. We define the struc- 
tural relaxation time, r a , as in ea. l|III.5|) . We find that, 
although the amplitudes of relaxation defer depending 
on k, the shear dependence of the relaxation time and 
its shear dependence is almost independent of k. The 
shear rate dependence of the relaxation time is given by 




FIG. 13: F(k,t)/S(k) for ka = 12. All symbols are as in 
FigEIl 



T a{i) ~ 7 ; consistent with the simulation results dis- 
cussed in the previous subsection. The four curves for a 
fixed shear rate but for different wave vector directions 
exhibit almost perfect isotropy, in qualitative agreement 
with the behavior observed in the simulations. Thus the 
apparent anisotropy of the vertex function in ea.l jll.43|l 
provides virtually no anisotropic scattering. For ka = 3, 
we see perfectly isotropic scattering. For ka = 6.7 and 
12.0, the curves show small anisotropy that is still con- 
sistent with the simulation results. Note that the dif- 
ferences shown in the plateau value for ka = 12.0 is an 
artifact due to our use of a square grid which produces 
an error in the radial distances depending on direction. 
The differences are noticeable but small in (b) and (d) 
for ka = 12.0. In these calculations, we have used the 
static structure factor at equilibrium and the distortion 
of the structure due to the shear was neglected. How- 
ever, as shown in Figures 3-6, the shape of S(k) becomes 
weakly anisotropic under shear. In order to see if the 
anisotropy of the structures affects the dynamics, we have 
implemented the same MCT calculation using S(k) and 
nc(fc) with an anisotropic sinusoidal modulation mimick- 
ing those observed in the simulations. We found that, as 
long as modulations are small, no qualitative change was 
observed and the dynamics was still isotropic. 

It is surprising that, although the perturbation is 
highly anisotropic, the dynamics of fluctuations are al- 
most isotropic. The reason for the isotropic nature of 
fluctuations may be understood as follows: The shear 
flow perturbs and randomizes the phase of coupling be- 
tween different modes. This perturbation dissipates the 
cage that transiently immobilizes particles. Mathemati- 
cally, this is reflected through the time dependence of the 
vertex. This "phase randomization" occurs irrespective 
of the direction of the wave vector, which results in the 
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essentially isotropic behavior of relaxation. This mecha- 
nism is different from that of many complex fluids and 
of dynamics near a critical point under shear, in which 
anisotropic distortion of the fluctuations at small wave 
vectors by shear plays an essential role J37j. 



D. Viscosity 

The shear-dependent viscosity 77(7) is evaluated by 
modifying the mode-coupling expression for the viscos- 
ity near equilibrium [3]. This can be done by following 
the same procedure explained in Section |H] for the to- 
tal momentum both from solvent molecules and colloids. 
The equation for the total momentum is given by the 
Navier-Stokes equation. The nonlinear term in the equa- 
tion comes from the osmotic pressure of the suspension 
particles. Neglecting the coupling of the density field of 
colloids with that of the solvent molecules, one obtain the 
mode-coupling expression for the viscosity. In the pres- 
ence of the shear, the same modification as in ea. (|11.37() 
is necessary and the wave vectors should be replaced by 
their time-dependent counterparts. The final result is 
given by 



1 

W 
dk 



k x k x {t) 



dS(k) dS(k(t)) ^ 



(2tt) 2 S 2 (k)S 2 (k(t)) dk y dk y {t) 



(IV.8) 

where 770 is the viscosity of the solvent alone. The integral 
over k can be implemented for the set of F(k, t) evaluated 
using ea. I|ll. 421) . In Fig. 14, we have plotted the shear 
dependence of the reduced viscosity defined by 



= vii) - vo 



(IV.9) 



for various densities around <f> c . The strong non- 
Newtonian behavior is observed at high shear rate and 
large densities, which is again in qualitative agreement 
with the simulation results for liquids reported in Ref. [|| . 
The shear thinning exponent extracted from the data be- 
tween 1CT 10 < Pc < 1 is i] R (j) 7 -!y with v ~ 0.99. 
This is in agreement with the exponent estimated from 
the simulation results and the structural relaxation time, 
T a (j), observed in the simulation and in the theory dis- 
cussed in the previous subsections. For larger shear rates, 
Pe > 1, the exponent becomes smaller, which is again 
consistent with computer simulation 0. However, one 
should not trust this reasoning for Pe >1. As discussed 
in Subsection IIV Al it is expected that the distortion of 
structure c(k) and S(k) by shear becomes important and 
one should take these effects into account in the theory. 
In this regime, however, the glassy structure has already 
been destroyed and the system become more like a "liq- 
uid" . In the liquid state, the shear-thinning exponent is 
always expected to be smaller than 1 and range between 
0.5 and 0.80, 0, S IH . 




Pe 



FIG. 14: The reduced viscosity, ea. l|IV.9(l . versus the shear 
rate for various densities. Here Pe = 7a 2 /Do- From top to 
bottom; <f> = 0.766549, 0.76650, 0.766453, 0.76640, 0.76600, 
and 0.75600. The highest density is 4 x 10~ 5 % larger than 



Slightly above </> c , plastic behavior is observed, which 
implies the presence of the yield stress. Note that this is 
an artifact of the theory because (f> c predicted from the 
mode-coupling theory is much lower than the real glass 
transition density. Below the real glass transition density 
<f> g , it is expected that shear thinning behavior similar to 
that predicted by MCT for (j> < <f> c will result. 



V. CONCLUSIONS 

In this paper, we presented the derivation of the mode- 
coupling equation for the realistic supercooled fluids un- 
der shear and compared the results with the molecular 
dynamic simulation. Our starting point is fluctuating 
hydrodynamics extended to the molecular length scales. 
A simple closed equation for the intermediate scatter- 
ing function was derived. We applied the theory to a 
two-dimensional colloidal suspension with hard-core in- 
teractions. The numerical analysis of the equation re- 
vealed very good agreement with simulation results for 
a related system: a binary liquid interacting with a soft- 
core potential. Theory and simulation showed common 
features such as (i) drastic reduction of relaxation times 
and the viscosity and (ii) nearly isotropic relaxation ir- 
respective of the direction of the shear flow. The fact 
that the dynamics is almost isotropic supports the valid- 
ity of the schematic models proposed so far, in which the 
anisotropic nature of the nonequilibrium states was not 
explicitly considered [pi Il7| ■ 

The mode-coupling theory developed in this paper is 
far from complete. The most crucial approximation is 
the use of the 1st FDT, which was employed when we 
close the equation for dynamical correlators. It is al- 
ready known that the 1st FDT is violated for supercooled 
systems under shear as well as aging systems pj. With- 
out the 1st FDT, one has to solve simultaneously the set 
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of mode-coupling equations for the propagator and cor- 
relation function, which couple each other through the 
memory kernels. Research in this direction is under way. 
A simple argument, however, should suffice to explain 
the success of the present theory with regard to shear 
thinning behavior. Consider a system that evolves out of 
equilibrium with one effective temperature T e g. Crudely, 
we can use the effective steady state version of FDT vi- 
olation in the form 

0(t) dC ZJ (t) 



Xij (t) 



fc R T f 



B J eff 



dt 



(V.l) 



to eliminate the response function Xij{t) m favor of the 
correlation function (t) from the set of mode-coupling 
equations. We thus expect for the behavior of the struc- 
tural relaxation time, T a (j), will be unaffected as long as 
ea. (|V.lfl holds. Another important approximation was 
to neglect the small distortion of the structure (c(k) and 
S(k)) due to shear. The simulation supports that the dis- 
tortion is negligibly small at the low shear rates (small 
Peclet number). The construction of the equation for 
equal-time correlation functions such as the structure 
factor might be more subtle and should be considered 
in future. Simulation results show that the structure is 
distorted in a noticeable and anisotropic manner at high 
shear rates; the peak height of the first peak of S(k) was 
lowered as much as 10 percent at 7 = 10 . ft is inter- 
esting that the dynamics is still isotropic and qualitative 
behavior is not affected even for such high shear rates. 

In the absence of shear, MCT is known to break down 
at densities well below the real glass-transition density, 
where MCT incorrectly predicts a fictitious non-ergodic 
transition. Beyond this MCT cross-over density, acti- 
vated hopping between the local minima of the free en- 
ergy surface is expected to dominate the dynamics of the 
system. Sollich et aLp^. |47| have analyzed a schematic 
model for hopping processes and explained similar shear 
thinning behavior. Lacks [48j has also rationalized shear 
thinning behavior in terms of changes of the free energy 
barriers due to shear. It is interesting that the totally dif- 
ferent picture given by MCT leads to some qualitatively 
similar conclusions. The analysis of the violation of the 
1st FDT and effective temperatures has the possibility of 
clarifying the difference between these distinct pictures. 
Barrat et al. have shown by simulation that the 1st FDT 
is violated in a sheared liquid and observed non-trivial ef- 
fective temperatures 0- This is consistent with the con- 
clusion of Berthier et al.\l4fy for the non-equilibrium p- 
spin model. On the other hand, the trap model predicts 
multiple effective temperatures |47| . Future effort will be 
directed towards extracting effective temperatures from 
our fluctuating hydrodynamic approach. 

In Section II we have mentioned that there are sub- 
tle problems in constructing an equation for the inter- 
mediate scattering function alone: The "correct" MCT 
equation can be derived if you start from the nonlinear 
Langevin equations for the density and momentum fields, 
taking the overdamped limit at the end. But difficulties 
arises if the overdamped limit is taken for the Langevin 



equation at the beginning. These problems exist even for 
systems at equilibrium and are related to generic prob- 
lems that exist in loop expansion methods. These sub- 
tleties were already recognized in the mid-70's[34| and 
recently pointed out in the context of glassy systems by 
Schmitz et al. 0] . A detailed study of these and issues re- 
lated to the field-theoretic approach to out-of-equilibrium 
glassy liquids will be discussed in a future publication|32j. 
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APPENDIX A: STATIC STRUCTURE FOR THE 
TWO-DIMENSIONAL FLUIDS WITH THE 
HARD-CORE INTERACTION 

Baus et al. have derived an approximated analytic ex- 
pression of the direct correlation function c(k) for the 
d-dimensional hard sphere fluids |43|. For the two dimen- 
sion system, it is given by 



aka/2 



+ v I] dx (1 - 2:2)1/2 (n^ " iax)2 ^S^-) 

(A.l) 

where J n (x) is the Bessel function of the first kind, <f> — 
7r<7 2 /?o/4 is the packing fraction, po is the number density, 
and a is the diameter of spheres, a is a parameter which 
is determined by solving 



a 2 (a 2 - 4) arcsin Q j - (a 2 + 2) V 'a 2 - 1 



V 



(A.2) 



Z((f>) = p/pok^T is the compressibility factor which is 
expanded by a rescaled virial series as 



z(4>) = (1 



(A.3) 



where c n is a coefficient which is related to the virial 
coefficients via a recurrence relation. We truncated the 
series at n = 6. The coefficients c„ are given by c\ = 0, 
c 2 = 0.1280, c 3 = 0.0018, c 4 = -0.0507, c 5 = -0.0533, 
and c 6 = -0.0410. 
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